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ABSTRACT 

Mounting discoveries of extrasolar planets orbiting post-main sequence stars mo- 
tivate studies aimed at understanding the fate of these planets. In the traditional 
"adiabatic" approximation, a secondary's eccentricity remains constant during 
stellar mass loss. Here, we remove this approximation, investigate the full two- 
body point-mass problem with isotropic mass loss, and illustrate the resulting 
dynamical evolution. The magnitude and duration of a star's mass loss com- 
bined with a secondary's initial orbital characteristics might provoke ejection, 
modest eccentricity pumping, or even circularisation of the orbit. We conclude 
that Oort clouds and wide-separation planets may be dynamically ejected from 
IMq — 7Mq parent stars during AGB evolution. The vast majority of plane- 
tary material which survives a supernova from a 7Mq — 20Mq progenitor will 
be dynamically ejected from the system, placing limits on the existence of first- 
generation pulsar planets. Planets around > 2OM0 black hole progenitors may 
easily survive or readily be ejected depending on the core collapse and super- 
wind models applied. Material ejected during stellar evolution might contribute 
significantly to the free-floating planetary population. 

Key words: planet-star interactions, planets and satellites: dynamical evo- 
lution and stability, stars: evolution, stars: AGB and post-AGB, Oort Cloud, 
supernovae: general 



o 



X 



1 INTRODUCTION 

Understanding the formation and subsequent dynami- 
cal evolution of exoplanets has been a motivational hall- 
mark for many observational and theoretical investiga- 
tions. However, extrasolar planets continue to be discov- 
ered in surprising and exotic environments, and questions 
about the endstate of exoplanets are becoming increas- 
ingly relevant. Few studies so far have modeled these 
systems, which often feature evolved and variable par- 
ent stars. The rich dynamics therein fundamentally differ 
from studies of planets around main sequence stars. 

Examples of exoplanets which do not orbit main 
sequence stars are growing. The first confirmed ex- 
trasolar planets were discovered around a neutron 
star: specifically, the m il lisecond pulsar PSR1257-I-12 
(jWolszczan fc Fraill [l993 : IWolszczanI \im4 ). The mini- 
mum masses of these three planets continue to be among 
the lowest known to date, and two of these planets reso- 
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nantly interact. ISigurdsson et al.l l|2003h later discovered 
another pulsar planet, around the binary radio millisec- 
ond pulsar PSR B1620-26. Exoplanets are also thought 
to orbit white dwarfs and stars with wh ite dwarf com- 
panions. In the first category, GD 66 (IMullallv et al.l 
l2008l . |2009| ). GD 356 (|Wickramasinghe et al.ll2010l ) and 
Gliese 3483 (Matt Burleigh, private communication) are 
planet-hosting stars. In the second cat egory, examples are 
thought to include Gl 86 = HP 13445 JOueloz et al.[l2000l : 
iMugrauer fc Neuhau sei'12005 1: iLagrange et al.ll200g) . HD 
27442 llButler et al.l 2001: C hauvin et all 120061) .and HP 
147513 (JMavoret al.ll2004l : IPesidera fc BarbierillJOOTl ). 



Additionally, planets have been discovered orbiting 
stars that have turned o ff of the main seque nce but are 
not yet stellar remnants. ISilvotti et al.l (|2007|) discovered 
a giant planet orbiting the e xtrem e horizontal branch star 
V 391 Pegasi, iGeier et al.l l|2009f ) fo und a planet aro und 
the hot subdwarf star HP 149382, iLee et al.l (|2009l) re- 
ported circumbinary pl anets to the sdB -| -M ec lipsing sys- 
tem HW Virginis, and ISetiawan et al.l l|2010l ) suggested 
that the planet orbiting the red horizontal branch star 
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HIP 13044b might be of extragalactic origin. Cataclysmic 
variables are another class of systems which might har- 
bor planets, and recently, planets around the cataclysmic 
variab les QS Vir (Qian et al. 2010a), DP Leo (Qian ct al. 
l2010bl ') and HU Aqr (JQian et al.1 l201ll ) have been re- 
ported. 

Prospects for discovering additi o nal planets orbit- 
ing white dwarfs iJDrake et al.1 l2010l : iFaedi et al]|201lh 
and extreme horizontal branch stars (Bear & SokeiJ 
I2OIII ) are promising, and observati onal campaigns to 
do so have already been initiated ( Hoean et al.l I2OO9I : 
iBenatti et aDl2010l : ISchuh et al.llioiol '). The Kevler mis- 
sion can detect even sm aller bodies around white dwarfs 
(JDi Stefano et al.ll2010l l. 

Theoretical investigations regarding the evolution 
of planets around post-main sequence stars have fo- 
cused primarily on planet engulfment and interaction 
with the expanding stellar enve lope, both for exoplanet s 
and specifi c ally f or th e Earth. Villaver fc Livid (|2007|), 
iMassarottil (|2008l ) and IVillaver fc Livid (|2009l ) use par- 
ticular stellar evolutionary tracks to determine ranges 
of semimajor axes at which planets are likely to be en- 
gulfed. In this regime, tidal modelling has a significant 
effect on the s ubsequent orbi tal evolution. However, as 
summarized bv iHansenl ((20101), the nature of tidal dissi- 
pation is poorly understood and continues to yield dif- 
ferent results depending on the model and assumptions 
used. For this reason , the f at e of the Earth is uncer- 
tain. ISackmann et al.1 (119931). iRvbicki fc Denli l|200ll '). 



ISchroder fc Connon SmithI (|2008l ) and llorid (|2O10l ) all ( 
plore the fate of the Earth in light of the Sun's post 
main- sequence mass loss, w i th dif fering results. Alterna- 
tively, iDebes fc Siexirdssonl (|2002l ) focus on the stability 
of multi-planet systems and link stellar mass loss to in- 
stability timescales. By doing so, they demonstrate how 
multiple planets beyond the reach of the star's expanding 
envelope might become unstable. 

In this study, we consider just a single planet, or 
smaller body. We perform a detailed analysis of the vari- 
able mass two-body problem and apply the results to a 
wide range of star-planet fates that encompass all stel- 
lar masses < ISOA^q. We focus on how stellar mass loss 
affects the eccentricity of a planet or planetary material, 
a link often ignored in previous studies. As a result, we 
show that planetary material can be ejected from a sys- 
tem based on mass loss alone. We then quantify for what 
combination of parameters we can expect this behavior. 

We start, in Section 2, by reviewing the history of the 
variable mass two-body problem and the corresponding 
equations of motion. We then analyze the orbital evo- 
lution in different mass loss regimes, determine where 
and when the traditionally-used adiabatic approximation 
holds, and estimate when the planets would become un- 
stable. In Section 3, we apply the theory to stars of all 
masses up to ISOM© in order to pinpoint realistic sys- 
tems which would yield instability. We treat five different 
mass regimes in separate subsections. We then discuss the 
caveats, implications and potential extensions in Section 
4, and conclude in Section 5. 



2 THE GENERAL TWO-BODY MASS-LOSS 
PROBLEM 

2.1 Overview 

Mass loss in the two-body problem, where both bodies 
are considered to be point masses , has boon studied for 
over a century (e.g. lGvldenlll884l : [Mcstschcrsky 1893). 
This situation is sometimes referred to as the "Gylden- 
Mestschersky" problem, even though this particular case 
refers to both variable mass rates having the same func- 
tional dependence. The more general problem takes many 
forms, or special cases, which are nicely outlined by Ta- 
ble 1 of lRazbitnaval (|l98a l. One well-known form results 
from the application of this general th eory to binary stel- 
lar systems, a formalism pioneered bv lJeand (|l924l ). The 
mass loss prescription which bears his name, M = —kM-' , 
where M is mass and «; and j axe constants, has been an- 
alytically and numerically treated in many subsequent 
studies. However, specific applications of mass loss to 
planetary systems have received little treatment. 

Soon after the ad vent of co mputer-based numeri- 
cal integrations, iHadiidemetrioul (|l963l . [l966a|[b|) revis- 
ited and reform u lated the problem in important ways. 
IHadiidemetrioul (| 19631 ') highlighted the subtlety with 
which mass loss must be treated in order to retain phys- 
ical interpretations of the evolution of orbital elements. 
He modeled mass loss as an additional acceleration that 
is a function of a time- and mass-dependent velocity, and 
showed that for any isotropic mass loss prescription, a 
planet's angular momentum h satisfies: 



h — constant — \/Gfia (1 



(1) 



where a refers to the semimajor axis, e to the eccen- 
tricity, and n = Mi, + M. The subscript "•" refers to 
the star and the variables without subscripts refer to the 
(lower-mass) secondary in the two-body system, which 
can be thought of as either a planet or particle; we will 
use the term "planet." Despite the conservation of angu- 
lar momentum, no such conservation clai m could be made 
about the total energy of the system, rl IHadiidemetrioul 
(|l966bl ) then significantly discovered that amidst great 
mass loss, such as in a supernova, the eccentricity of 
the secondary may increase, and eventually lead to ejec- 
tion from the system. That finding is the foundation for 
this work. A s ubsequent series of papers (jVerhulstll 19691 : 
IVerhulst fc"Ec khaus 1970; Verhulst 1971) provided an ex- 
pansion of and comp arison with Hadjidemetriou's results. 
lAlcock et al.l (|l98d ) then approached the ejection possi- 
bilities from a different perspective by considering the 
effect of vigorous mass loss of white dwarf progenitors on 
a comet. Later, iParriott fc Alcockl (|l998l ) demonstrated 
how the asymmetric mass loss case yields a greater frac- 
tion of cometary ejections. 

Despite the wide body of work on mass loss in the 



^ Because a system with isotropic mass loss will maintain its 
rotational symmetry, according to Noether's Theorem, the an- 
gular momentum will be conserved. Because the same system 
does not exhibit time invariance, the energy of the system is 
not guaranteed to be conserved. 
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two-body problerrjj, most studies continued to concen- 
trate on binary stars. IPebes fc Sieurdssonl ()2002l ) helped 
break this trend by analyzing the planetary case through 
the modelling of multiple planets orbiting a single star. 
They assumed the planets had equal masses and initially 
circular orbits, and studied their motion in the "adia- 
batic" approximation. This approximation holds when 
the mass loss timescale is much greater than a plane- 
tary orbital timescale. In this approximation, the planet's 
eccentricity is thought to remain nearly constant, and 
hence, from Eq. ([l|, 



da\ 

'^Vadiabatic 



a dfi 
fi dt 



(2) 



However, in the general planetary case, the angular 
momentum is a function of eccentricity, which is generally 
not constrained to be fixed. Other complicating factors 
are: i) because planetary orbits which are changing due to 
stellar mass loss are not closed, averaged orbital element 
expressions can be misleadi ng and coun ter-intuitive, al- 
though technically correct (|lorid l20ld ). ii) in a single 
phase of stellar evolution, mass loss is typically noncon- 
stant (although monotonic) and may not be isotropic, iii) 
stellar mass evolution typically involves multiple phases 
of mass loss on timescales which can vary by orders of 
magnitude, and iv) several additional forces due to stel- 
lar evolution, such as tides and dynamical friction from 
the expanding envelope, might be necessary to model in 
order to describe the correct orbital evolution. 

Here, we do not place restrictions a planet's semima- 
jor axis, eccentricity or orbital angles, but do take mea- 
sures to focus our results. We treat mass loss as isotropic. 
Tidal effects are unimportant in the regimes we consider 
here, and so we can safely neglect those. To foster intu- 
ition for the mass loss problem, and to obtain tractable 
results, our analytics assume a constant mass loss rate 
throughout. However, some of our analytical results are 
completely independent of the mass loss rate assumed. 
The parameters for the example cases used in this Section 
were selected to best demonstrate different aspects of the 
motion of this general two-body problem with mass loss; 
more realistic cases are presented in Section 3. There, we 
apply the theory presented here to just a single phase of 
stellar evolution, but do consider almost the entire phase 
space of stellar mass. 



2.2 Statement of Equations 

Although the equations of motion in terms of orbital el- 
ements for the variable-mass two-body problem can be 
derived from fi r st pr i nciples, only a few authors ( e.g. 
Hadiidemetrloul Il963l : IVerhulstI 1 19691 : iDeprid Il983l : M 
20081 ) have stated them in full without averaging or ap- 
proximation: 



^ iRahoma et al.l 1 120091 ) provides a detailed summary of addi - 
tional results from past papers, and lPlastino fc Muzziol 1119921 ) 
summarizes the "use and abuse" of using a force to model mass 
loss. 
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dt 
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dn 

dt 
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~dt " 
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dF 



sin / 1 dfj, 
e fi dt 
n{l + ecos/) 



(1 



„2^3/2 



(3) 

(4) 

(5) 
(6) 

(7) 



where i is the inclination, Q, is the longitude of ascending 
node, uj is the longitude of pericenter, uj is the argument 
of pericenter and / is the true anomaly. Equations fi]), ((3| 
and (U) are self-consistent and may be derived from one 
another with help from the vis-viva equation. The time 
derivative of position in terms of orbital elements and the 
statement of the conservation of angular momentum in 
polar coordinates give Eqs. ((6]) and ^. 

These equations may also be derived from more gen- 
eral considerations. Gauge theory is a basis from which 
one may obtain sets of equations such as Lagrange's plan- 
etary equations and Gauss' Planetary Equations by defin- 
ing just a single perturbative acceleration to the clas- 
sic two-body problem, and a gauge velocity. The for- 
mulation of the theory with regard to planetary dy- 
namics as well as e xtensive descriptions can be foun d 
in Efroi msky fc Gol drcich (2003, 2004), Gurfil (|2004|) 
Efr oimskvn2005al). lEfroimskyl (l2005bl |2006|). iGurfill 
(200^) and iGurfilfc BelvaninT (|2008l ). iHadiidemetrioul 
(1963J) showed that the sum of the isotropic mass varia- 
tion of both bodies is equivalent to a perturbative force 
with an acceleration of AA = — (l/2)(d/i/df)(l//i)v, 
where v is velocity. This acceleration yields Eqs. ©-0 
directly for a zero gauge. 

Every variable in Eqs. ©-(O is considered to be 
a function of time. The mean motion, n, is equal to 
G^'^fi^'^ /a"^'^, where G is treated as the standard grav- 
itational constant. Although we use /i throughout this 
work to emphasize how the motion is affected by the sum 
of the mass loss (or gained) by both bodies, the value of 
the planetary mass and how it changes with time has a 
negligible effect on the results for A/* ^ M. For a IMq 
star, if one assumes a planetary mass of ~10 Jupiter 
masses, which is on the order of the theoretical upper 
bound, then M/M^, ~ 1%. 

The planet's true longitude, 0, varies according to: 



d9 n{l + ecos/) 
rfi" (l-e2)3/2 ' 



(8) 



which is not explicitly dependent on the mass loss rate 
and hence is equivalent to the case of no mass loss. This 
equation demonstrates that from the point of view from 
a fixed reference direction, the secondary will continue to 
circulate around a star that is losing mass as long as the 
secondary remains bound. 

For completeness, we consider the evolution of other 
traditionally-used orbital parameters. The planet's eccen- 
tric anomaly, E, will vary according to: 
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dE n(l + ecos/) 
'dt " l-e2 



sin / 1 dfi 



e\/l — e2 /i di 



(9) 



Note that the right-hand sides of Eqs. ((3])-(|9]) may be 
expressed in terms of the eccentric anomaly instead of 
the true anomaly. The planet's mean motion will vary 
according to: 

dn n (2 + e^ + 3e cos /) 1 dp 
dt 1 — e^ n dt 

The planet's mean anomaly, 11, can be expressed as an 
explicit function of time by use of the "time of pericen- 
ter", r: 



(10) 



dn , , (2 



l-e2 
or, through Kepler's Equation, as 

dn 



3e cos /) 1 dp 
/i dt 



(11) 



= n + 



Vl — e^sin/ (l + e^ + ecos/) 1 dfi 



--^^ (12) 



dt ' e(l + ecosf) fidt' 

which is explicitly independent of time. Finally, the mean 
longitude. A, changes with time according to 



dX dn dzu 

dt dt dt 

duj 

n — 

dt 



Vl-e2(l + e^ + ecos/) 
1 + e cos / 



-1 



(13) 



Throughout this paper, we denote initial values with the 
subscript "0". 

2.3 Parametrizing Mass Loss 

Suppose the mass loss rate is constant and equal to —a, 
such that Q > 0. Then fi = G (po — <^t), and 



fi dt J \ a J 



(14) 



We can better quantify adiabaticity and various 
regimes of motion due to mass loss by defining a dimen- 
sionless "mass loss index", ^: 



* = — 



a 
np 

1 

2^ 



t) VIAU/ 



M 



IMq/jtJ VlAU/ \IMqJ ■ ^^^' 

This parameter provides a scaled ratio of the orbital 
period to the mass loss timescale. The initial value of the 
index as ^q. nence, the time evolution of 'I' is governed 
by: 

1 + e cos / ^ 1 d/i 
fi dt 



— - = -3* 

dt 



1 



(16) 



When $ ^ 1, a system can be considered "adiabatic" 
the case we treat first. 



2.4 "Adiabatic" Regime Evolution 

2.4-1 Adiabatic Eccentricity Evolution 

We begin analyzing the equations of motion by first con- 
sidering Eq. (|4]), because all nonzero equations of motion 
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Figure 1. Analytic approximation to the eccentricity evolu- 
tion in the adiabatic regime. Shown here is the difference in 
eccentricity of a planet evolving according to Eq. I I17I I com- 
pared to Eqs. ©-([Til. The planet is located at ao = 1 AU 
from a /iQ = lAfQ star which is losing mass at the rate of a 
rate of a = IQ-'^Mq/yt (*o Ri 1.6 X lO""). The differently 
coloured lines from the top of each crest moving downward 
correspond to eo = 0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8, and 
0.9 respectively. 



explicitly include e in some manner. Note importantly 
that the equation demonstrates that an initially circular 
planet will not remain on a circular orbit, and that the 
planet's eccentricity will undergo oscillations on orbital 
timescales when the parent star loses mass. 

We can solve Eq. (|3| by noting that in the adia- 
batic approximation (^ <C 1), the first term in Eq. (O is 
considered to be negligible compared to the second term 
{— dO/dt), because the first term is proportional to the 
mass loss rate. Further, /i is assumed to remain fixed over 
the course of one orbit. Hence, in this regime, Eq. Q may 
be integrated directly over the true anomaly, with the re- 
sult: 



eadiabatic = Co + *0 - 



l-eg)^sin/ 
1 — eo cos / 



(17) 



According to Eq. (|17p . after each orbit the eccentric- 
ity will return to its initial value. During the orbit, the 
amplitude of (eadiabatic - eo) is *o(l - Gq) cx a. Thus, 
assuming a current value of Qq « 10~^^/yr, the Earth's 
eccentricity is raised by about 10^^^ each year due to the 
Sun's mass loss. 

Fig. [T] demonstrates the accuracy of Eq. p7[) when 
compared with the evolution from the full equations of 
motion (Eqs.[3][7|) for a, Oo = 1 AU planet orbiting a po = 
IMq star which is losing mass at the rate of lO~^M0/yr 
(8 orders of magnitude greater than a©). The agreement 
is excellent over the course of a single orbit. Over time, 
the approximation gradually worsens, as the evolution of 
^pQ is not taken into account in Eq. (|17|) . 
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Figure 2. The adiabatic regime. The position in space {left panel) and the evolution of the true anomaly (right panel) of a 
planet (or belt particle) that is being pushed outward due to stellar mass loss. The colours on the curves indicate evolution at 
the same points in time, and the vertical lines of true anomaly indicate circulations of the angle. The star has an initial mass of 
^0 = ^* = 1^0 S'lid is losing mass at the rate of a = 1 X 10~^MQ/yr. The planet begins on a highly eccentric orbit (eo = 0.9) 
at ao = 100 AU, with /o = 0°. Hence, '^o = 0.0016. Notice that as the planet moves outward and its mean motion decreases, the 
circulation period of the true anomaly decreases as well. 



2.4-2 Adiabatic Semimajor Axis Evolution 

We now consider the semimajor axis evolution from Eq. 
(131). Note from the equation that for any mass loss, the 
semimajor axis can never decrease. 

In the adiabatic regime, the semimajor axis is tra- 
ditionally evolved according to Eq. ^. Note, however, 
that the equation does not follow from Eq. ((3| if e 7^ 0. 
Yet, when the semimajor axis is averaged over one or- 
bital period, the eccentricity terms vanish and Eq. is 
recovered. The solution of this equation is: 



^adiabatic 



■ ao 



at 
MO 



(18) 



Therefore, an adiabatically evolving planet will, for 
example, double its orbital separation if its parent IMq 
star constantly loses mass at the rate of a = 5 x 
W'^Mq/yt over 100 Myr. In a different example, a 2Mq 
star is expected to lose at most « 70% of its initial mass. 
Therefore, if all this mass is lost adiabatically, then orbit- 
ing planets can expect to increase their semimajor axis 
by at most a factor of ~ 3.3. 



2.4-3 Adiabatic Orbital Angle Evolution 

Turning to other orbital parameters, the longitude of 
pericenter is a typically secular feature of multi-planet 
extrasolar systems. Its variational timescale is often on 
the order of thousands of orbits. During stellar evolu- 
tion, however, Eq. ([6]) demonstrates that the variation 
in a planet's longitude of pericenter is quick (on orbital 
timescales), and changes sign over each orbital period. To 
be consistent with the adiabatic approximation, in which 
duj/dt Ri in Eq. 0, then 



Because dvj/dt is assumed to be zero, the value of 
/adiabatic foUows the sauic evolutiou as / would in the 
two-body problem with no mass loss. 

We can obtain an adiabat for n from Eq. (|10|) under 
the same assumptions that were used to derive Eqs. ((2| 
and (fUl): 



(20) 



^adiabatic — J^O 1 

MO/ 



'^adiabatic 



■ ^0- 



(19) 



Thus, in the adiabatic approximation, the mean motion 
is a monotonically decreasing function. In the same ex- 
ample system from Section 2.4.2, with /lo = IMq, and 
Of = 5 X 10~^MQ/yr, after t = 100 Myr the planet's mean 
motion would decrease by a factor of 4. This result is ex- 
pected from Kepler's 3rd law with a halved stellar mass a 
doubled semimajor axis. For 2Mq stars, the final Keple- 
rian period of a planet when mass loss has ceased would 
be enhanced from its initial period by a factor of at most 
^ 11. 



2-4-4 Adiabatic Evolution in Space 

In space, adiabatic evolution corresponds to a planet or- 
biting in an outward spiral pattern. Figure[2]displays such 
an orbit (for ^o ~ 0.0016), which is not closed. After each 
cycle of true longitude, the eccentricity does return to its 
initial osculating value. The semimajor axis is seen to in- 
crease by as much as 10% of ao per orbit. The increase in 
orbital period can be linked with the circulation timescale 
of/. 

A highly eccentric planet might make close passes 
to the star, close enough to be affected by tides and the 
stellar envelope. In order to determine if the planet is 
more or less likely to suffer these encounters from mass 
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loss, consider the evolution of the pericenter, q, of the 
planet: 



dq 
'di 



a(l-e)(l 



1 + e 



cos /) 1 d/i 
/i dt 



(21) 



Equation ((21} demonstrates that the pericenter mono- 
tonically increases with stellar mass loss. The left panel 
of Fig. [2] corroborates this relation. Therefore, if a planet 
"outruns" the star's expanding envelope, then one can 
neglect the envelope's influence on the planet. 



2.5 Regime Transition 

2. 5. 1 The Breaking of Adiabaticity 

Equation (|20p has important implications for the dynam- 
ical system, as mean motion is inversely proportional to 
the Keplerian period. Hence, as a star loses mass, and 
pushes a planet radially outward, the mean motion de- 
creases, and eventually the orbital period will be compa- 
rable to the mass loss timescale (*1/ ~ 1). More precisely, 
dd/dt (which is proportional to the mean motion) will 
eventually become equal to {—dzu/dt) (which is propor- 
tional to the mass loss timescale). At this bifurcation 
point in the dynamics, Eq. (O demonstrates that the true 
anomaly becomes momentarily stationary. At this point, 
one can claim that adiabaticity is broken. 

Note that in the adiabatic regime, / circulates. At 
and beyond the bifurcation point, df /dt instead begins 
to librate. The effect of a librating / on da/dt and 
de/dt is pronounced, quick and runaway. The eccentric- 
ity and semimajor axes evolution undergo a qualitative 
change, and the rate of increase in the latter is pro- 
nounced. Therefore, we denote this regime as the "run- 
away" regime. We wish to investigate this regime tran- 
sition, and do so first qualitatively through Figs. [3] and 
m These figures model a IMq star which is losing mass 
at a relatively high rate of q = 5 x lO~^M0/yr over 
1.5 X 10* yr. After this amount of time, the star will have 
lost 75% of its initial mass. These values are chosen for 
demonstration purposes, as typical IMq stars will lose 
~ 35% — 62% of their mass en route to becoming a white 
dwarf. We model more realistic systems in Section 3. 

Figure [3] illustrates the approach to and onset of 
adiabaticity breaking, with 4'o = 0.023 (left panel) and 
^0 = 0.25 (right panel). In the first case, the eccentricity 
ceases to remain approximately constant and can start 
to oscillate on the order of a tenth. In the second case, 
where ^ quickly assumes values on the order of unity, 
planets evolve in the runaway regime and may achieve 
hyperbolic orbits. 

Figure U showcases the semimajor axis evolution for 
the same systems in Fig. [3] In the left panel of Fig. |4l 
the curves of initial eccentricity break away from the adi- 
abat, increasing at a steeper rate than the adiabat. In 
the right panel, at t = 0, the systems are just beyond 
the adiabat and are "running away" from the star. For 
a constant mass loss that turns on and off nearly instan- 
taneously (such as in a supernova), planets might not 
ever evolve adiabatically, and begin their life in the run- 
away regime. Note that unlike eccentricity, the semimajor 
axis is always increasing, even when oscillating about the 



adiabat. In the runaway regime, the departure from the 
adiabat is drastic; the right panel shows that the planet 
may increase its semimajor axis by many factors before 
achieving a hyperbolic orbit. 

The resulting increase in a will cause the mean mo- 
tion term in Eq. to decrease further. In the limiting 
case where duj/dt ^ n, the libration amplitude will ap- 
proach zero, and the true anomaly will become nearly 
stationary. As a result, the cos / term in Eq. Q becomes 
constant, and de/dt becomes linear in e, causing a posi- 
tive feedback effect that is characteristic of the runaway 
regime. 

An orbit that is transitioning out of adiabaticity will 
not complete its final orbit around the star, as the true 
anomaly is no longer circulating. Figure [5] illustrates the 
resulting motion in space. The system selected is the 
highest eccentricity (eo = 0.9) curve from the left pan- 
els of Figs. [3] and H (*o = 0.023). The system reaches 
the bifurcation point at ~ 1.2 x 10* yr, as one can read 
off from the right panel. Now we explore this bifurcation 
point analytically. 



2.5.2 Characterising the Bifurcation Point 

The bifurcation point, as we defined in the last subsec- 
tion, is the first moment when dO/dt — —dzu/dt. At this 
moment, from Eq. (O, df /dt = 0, and 



*bif 



ebif (1 + ebif cos /bif )^ 
sin/Mf(l-egif)^/^ 



(22) 



For the majority of possible values of eccentricity 
and true anomaly, ^^bif ~ 0.1 — 1.0. There are an infinite 
number of triples (*I'bif, ebif, /bit) that satisfy Eq. (|22[) . 
We cannot determine any of these three values from the 
initial conditions, although one may approximate ebif ~ 
eo from the adiabatic approximation. However, at this 
point in the planet's evolution, ebif might differ by over 
0.1 from eo. 

The functional form of Eq. H22p suggests that for a 
given ^Jbif and a given ebif, there might be more than one 
value of /bif which satisfies the equation. We now investi- 
gate this possible multiplicity further by considering the 
extremities of ^hu with respect to ebif and /bif. There 
are six values of /bif which satisfy d^&bif/rffbif = 0, five of 
which are unphysical. The one physical solution is: 



/bif.n 



i-ym 



sgif 



2ebif 



(23) 



where 90° sj /bif ^ 270° . Let the value of *bif at /bif = 
./bif, mill be denoted as *I'bif,fmin. Then, for a given 'I'bif and 
a given ebif, the number of values of /bif which satisfy Eq. 
are: 



values of /bif if *bif < *bif,fiiiiii (24) 

1 value of /bif if *bif = *bif,fmin (25) 

2 values of /bif if *bif > *bif,fmiii. (26) 

The maximum value of ^bif.min (obtained in the 
limit ebif — >■ 1) is 4/(3\/3). Therefore, for any given 
*I'bif > 4/(3\/3) ~ 0.77, there are two possible values of 
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Figure 3. The breaking of adiabaticity. Shown is the eccentricity evolution over 10* years of evolution of ao = 200 AU planets (or 
belt particles; left panel) and ag = 10'^ AU planets (or belt particles; right panel) orbiting a fiQ = M* = IMq star losing mass at 
a rate of q = 5 X lO~^M0/yr. The initial mass loss index for the systems in the left and right panels are respectively Sto x 0.023 
and '^0 ^ 0.25, values which are close to the transition point in the dynamics between the adiabatic and runaway regimes. Here, 
/o = 0°. The lines with an increasing dash length represent eo values of 0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, and 0.9, respectively. 
In the right panel, all planets are ejected from the system within 10* yr except the two planets with the lowest initial eccentricity. 
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Figure 4. The breaking of adiabaticity for the same two systems in Fig. [3] In the left panel, note how the eccentric planets 
oscillate about the adiabat until reaching the runaway regime. In the right panel, the planets begin at t = just off of the adiabat, 
and quickly settle into the runaway regime. 



/bif which satisfy Eq. (|22|l . Figure [6] demonstrates Eqs. 
p4l)- ((26)) graphicaUy by plotting *bif vs. /bif for 10 val- 
ues of fibif (0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, and 0.9). 
Large dots mark where ^bif = ^bif.fmin- As adiabatic sys- 
tems increase ^ and approach the bifurcation point, their 
evolution can be imagined as moving upwards on this plot 
while circulating almost parallel to the X-axis. Eventu- 
ally they will reach the bifurcation point, preferentially 

at ^bit.min- 

Note from Eq. ([22]) that ebif -> implies *bif -> 0, 
suggesting that planets with initially circular orbits can 
never be in the adiabatic regime. However, this is not 
true. If eo = 0, then from Eq. (flT)) . Cadiabatic = 'I'osin/. 
Inserting this expression into Eq. Q yields: 



4fj 

dt lt = 0,eo=0 



(l-F^^osin2/)' 
(l-*gsinV)'^' 



>0 



(27) 



for any nonzero value of / (even if /o — 0, then / at- 
tains a positive value immediately). After i = 0, df/dt 
will then continue to increase until the bifurcation point 
is reached. Therefore, initially circular planets may easily 
evolve adiabatically, which corroborates numerical simu- 
lations. 

Now we consider d^bif/debif = 0. There are 3 solu- 
tions, 2 of which are physical: 



ebif.cxtl 



i ( -3 cos /bit + J -I + I cos (2/bif) |;28) 



ebif,cxt2 = - I -3 cos /bif ^ y ^9 + o cos (2/bif) (29) 



where /crit < /bif < (360° - /ciit), such that 
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/crit 



180 



1 -1 

-COS 



160.53 



(30) 



Solutions At The Bifurcation Point 



J are: 


X) 


0° < /bif < ./cot 


(31) 


/cot < ,/bif < 180° 


(32) 


/bif ~ ./crit 


(33) 


fibif = Ebif.oxti and 


Cbif = 



This critical true anomaly value will be important for de- 
scribing motion in the runaway regime because it deter- 
mines where a quahtative change in the evolution occurs. 
For a given ^bif and a given /bif, the number of values 
of Cbif which satisfy Eq. (1221) 



1 value of Cbif if 

3 values of ebif if 

oo values of ebif if 

Limiting values of ^bif at 
ebif,cxt2 are 2/3 and 4/(3\/3). 

We can illustrate the multiplicity suggested by Eqs. 
([3l|-(|33l) with Fig. [71 Plotted in Fig.[7]are six curves cor- 
responding to / = 1°, 5°, 10° (short-dashed blue curves), 
and / = 179.0°, 179.5°, 179.9° (long-dashed red curves). 
We display these curves because they approximate e — *1/ 
evolution tracks beyond the bifurcation point in the / — >■ 
0° and / — >■ 180° cases. For these two values of /, 
df/dt « (from Eq.O, and hence / remains constant as 
a function of time. In this case, a system will move along 
one of these tracks. Because ^ is always increasing, for 
the blue short-dashed curves, this system will increase 
its eccentricity until ejection. However, for the red long- 
dashed curves, the planet's eccentricity might decrease 
until reaching a critical point, when the increase in vp 
will break the constant / approximation and force the 
system off the track. The critical points (circles; peaks) 
along these tracks which appear between 1/2 < e ^ l/\/2 
are given by: 



*bif 



6 cos /bif + e 



24x/3sin/bif 
\/5 - 3 cos (2/bif) - ecos /bif, 



(34) 



and for those critical points (triangles; troughs) which 
appear between l/\/2 ^ e < 1,: 




60 90 120 150 180 

/bif (deg) 



Figure 6. Values of (\I'bif i ^bif i /bif) a-t the bifurcation point, 
where lines with increasing dash length represent ebif values 
of 0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8, and 0.9, respectively. 
Corresponding colour dots represent the value of SPbif ,fmin for 
a given Ebif, when just one value of /bif satisfies Eq. ((22}. 
Adiabatic systems approaching the bifurcation point would be 
traveling upwards on this plot while circulating nearly parallel 
to the X axis. 



*bif 



6 cos /bif + e 



24\/3sin/bifj 
(8-6 cos^ /bif + e cos /bif) 
(5 - 3 cos (2/bif) -f ecos /bif)'/'' 



(35) 



where e = -^/IS cos (2/bif) — 14. The two critical curves 
are shown on Fig. [7] as thin black lines. The limiting val- 
ues represented by the three triangles, which are not dis- 
tinguishable by eye from one another on the plot, occur 
at ebif < 1. 

Linking the orbital parameters at the bifurcation 
point with the initial system orbital parameters is dif- 
ficult because although the bifurcation point is well- 
defined, the adiabatic approximations begin to break 



© 2011 RAS, MNRAS OOO.ITH231 



The Great Escape 9 



Solutions At The Bifurcation Point 
0.5 



lOV 

10^ ■ ^ • 

^ 10^/^. — z^- 

s ' ^ — 

u 

io-i> 



/f 



3V3~ 







0.2 



0.4 0.6 

fbif 



0.8 1.0 



Figure 7. Values of ('I'bifi ^bifi /bif) ^t the bifurcation 
point, where the blue-short dashed lines, starting from the 
top down, represent /bif = 1°,5°,10°, and the red long- 
dashed lines, starting from the bottom up, represent /bif = 
179.0°, 179.5°, 179.9°. The peak and trough critical points for 
the red long-dashed curves are marked with red circles and tri- 
angles, respectively. The black curves are given by Eqs. I)34| l 
and I I35I I. The blue and red curves, separated by /crit (Eq. 
I30| l exhibit qualitatively different behaviors. The blue short- 
dashed curves and the rising portions of the red long-dashed 
curves represent evolutionary tracks beyond the bifurcation 
point, demonstrating that for fait < / < 360° — /criti the 
planet's eccentricity will experience an initial decrease beyond 
the bifurcation point. 



down before the bifurcation point is reached (see, e.g., 
Fig. [3]). Nevertheless, we can analytically estimate the 
semimajor axis at the bifurcation point by using the semi- 
major axis adiabat. Doing so gives: 



fflbif 
ao 



%' 



1 2 

e^if (1 + ebif cos /bif) 3 
(l-egif)Msin/bif)3 



M. 



1 -i _i 



(36) 



Equation (|36|l contains qualitative physics useful for un- 
derstanding when the system reaches the bifurcation 
point. The dependence on the initial stellar mass, ini- 
tial semimajor axis and the mass loss rate determine how 
prone a star is to reaching the bifurcation point and eject- 
ing its planet. For two planets in separate systems with 
the same ao, the parent star whose physical parameters 
yield a smaller value of abif is more likely to cause ejec- 
tion. Also, wide orbit planets are more prone to be ejected 
than smaller orbit planets. Unfortunately, the term in 
square brackets is unknown and cannot be bound with- 
out some assumptions on ebif and /bif. To be consistent 
with using the adiabat, one can assume ebif = eo. How- 
ever, by the time the system has reached the bifurcation 
point, the eccentricity could have already varied away 
from its initial value by at least a tenth. The value of /bif 
is the cause of greater uncertainty. 



2.6 "Runaway" Regime Evolution 

2.6.1 Runaway True Anomaly Evolution 

The unknown value of /bif largely determines how the 
planet will evolve past the bifurcation point. If the mass 
loss is great and sudden enough, then the planet will by- 
pass the bifurcation point altogether and immediately 
start evolving in the runaway regime. In this case, the 
planet's /o value is crucial to its evolution. The divided 
phase space structure of Fig. [7] correctly suggests that 
systems can behave quantitatively differently depending 
on their true anomalies. 

Consider Figs. [8] and |9] which illustrate the eccen- 
tricity evolution as a function of initial true anomaly 
for 'I'o — 0.089 (approaching the bifurcation point) and 
'if — 7.96 (runaway regime). This dependence is more 
complex around the bifurcation point ^ ~ 0.1 — 1 than 
after it (^t > ^ibif). Note importantly that the divided 
phase space structure in Fig.Qmanifests itself strongly in 
Fig. E] (at /crit in the runaway regime, highlighted by the 
dotted blue box, when every planet's eccentricity must 
experience an initial decrease), but not in Fig. |8] (be- 
fore the bifurcation point). Additionally, in Fig. [H] for 
at least a third of all possible initial /o values, the first 
planet ejected has a eo value in between the extremes 
sampled of 0.01 and 0.9. Contrastingly, in Fig. ^ in ev- 
ery instance the first planet ejected has either eo = 0.01 
or eo — 0.9. Further, the eccentricity evolution is nearly 
symmetric about / — 180° in the runaway regime, a ten- 
dency not exhibited in Fig. [S] This helps to demonstrate 
how complex the evolution can be when the system is nei- 
ther robustly in the adiabatic or runaway regime. Planets 
which begin their post-main sequence life already in the 
strongly runaway regime ('I' ^ 1; Fig.|9j experience more 
predictable behavior. 

Fortunately, in some cases we can analytically ap- 
proximate the evolution of orbital parameters in this 
regime. If a planet begins to evolve att = 0inaVl'2>l 
system with a value of / that is close to either 0° or 180° , 
then / is guaranteed to librate with a small enough am- 
plitude so that / may be treated as a constant. Figures 
[7] and in] suggest that the resulting behavior in each of 
the two cases will differ qualitatively. The latter figure 
illustrates that for / — 180°, immediately after circular- 
isation, the eccentricity evolution starts increasing and 
continues to do so up until ejection. Before quantifying 
this behavior analytically, we first attempt to explain the 
physical mechanism at work: 

At / ~ 0°, e will increase until the planet is ejected. 
There is no alternative evolutionary track. At / ~ 180°, 
the eccentricity will decrease until e — >■ 0. In this limit, 
\dzu/dt\ becomes large, forcing df/dt ^ 0. The true 
anomaly will then quickly sample other values. At all 
other values except 0°, df/dt 7^ 0. When / eventually 
samples 0°, it becomes stuck on that evolutionary track. 



2.6.2 Runaway Eccentricity Evolution 

In the runaway regime, when / = 0° or / = 180°, Eqs. 
(|31) and @ may be solved directly and analytically. The 
eccentricity evolution is then given by: 
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Figure 8. How true anomaly affects eccentricity evolution on the approach to the runaway regime. Shown is the eccentricity 
of a planet at ao = 500 AU over 1.5 X 10^ yr orbiting the same star (/iq = ^Mq and a = 5 X lO~^M0/yr, so *o ~ 0.089) 
as in the left panels of Figs. [3] and 2] as a function of /q. The lines with increasing dash length represent eo values of 
0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8, and 0.9, respectively. Note the dramatic sensitivity the initial true anomaly may have on 
the eccentricity evolution, and that for /o = 120° — 240°, the particle or planet which is ejected first is one with an initial eccen- 
tricity that is neither the highest nor lowest sampled. The evolution is not symmetric about /o = 180° , and a few of the initially 
eccentric planets will become circularised. 
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Figure 9. How true anomaly affects eccentricity evolution in the runaway regime. Shown is the eccentricity of a planet at ao = 10* 
AU for the situation in Fig. |8] Here, however, Nto ~ 7.96 and the eccentricity evolution is nearly symmetric about /o = 180° (as 
hinted at by the /o = 170° and /o = 190° cases). Circular orbits are approached at /o = 180°. The blue dashed box highlights 
the case /o = /crit- For /crit < /o < 360° — /crit, every planet, regardless of eo, is predicted to experience an initial eccentricity 
decrease. The eccentricity will later increase if the mass loss continues for a long enough time (which is not the case, e.g., for 
/o = 180° and eo = 0.8). 
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(38) 



In the /o = 0° case, the eccentricity will increase until the 
planet is ejected; in the /o = 180° case, the eccentricity 
will decrease until the planet achieves a circular orbit. 
Hence, the amount of mass remaining in a star at the 
moment of ejection, Mout and at circulation, Mcirc, are: 



Mout 

Mo 

/^circ 

Mo 



1 + eo 



eo. 



(39) 



(40) 



Equation (|39|l demonstrates that for '^ > ^I'bif and 
/o ~ 0°, a planet will be ejected before half of the star's 
mass is lost. Also, planets with larger initial eccentricities 
would be the first to be ejected. Equation (|40p demon- 
strates that a planet of any eccentricity may be circu- 
larised, and that nearly initially circular planets are the 
most likely to do so first. These equations may also be 
expressed as tout = Mo (1 ^ eo) /(2a) and fcirc = fioeo/a. 

After circularisation, the planet's true anomaly 
quickly becomes 0°, as described in the last subsec- 
tion. Then, the eccentricity evolves according to a "post- 
circular" prescription: 



epost — circular — 



-i 



= ^(1- 
M 



eo)-l. 



(41) 



The total time taken for a planet to circularise and then 
be ejected from a system is ftot = mo (1 + eo) /(2a). After 
this time, the amount of material the star has depleted 



is: 

Mtot 
Mo 



eo 



(42) 



Any planet which circularises before becoming ejected 
therefore must have a parent star that loses at least half 
of its mass. Additionally, in order for the most eccentric 
planets to be ejected, they require the star to lose all 
of its mass. The implications are that no belt of objects 
that are uniformly distributed in both true anomaly and 
eccentricity can all be ejected due to mass loss: Regardless 
of the value of *&, the highest eccentricity bodies at / ~ 
180° must survive. 



2.6.3 Runaway Semimajor Axis Evolution 

The semimajor axis evolution in the /o = 0° , 180° run- 
away regime is: 



^runaway | / — 180° 



= ao (1 + eo) 



-1 at 



(l + ^^o)-^ 



ao (1 + eo) 



(1-eo) 



(44) 



respectively. As one might expect, for initially circular 
orbits in the runaway regime, the semimajor axis evolu- 
tion is the same for /o — 0° and /o = 180°. Also, in the 
circular limit, we can compare the semimajor axis evolu- 
tion with what it would have been in the adiabatic limit. 
For a given ao, flrunaway/aadiabatic = (2m/mo — 1)"^, which 
holds until /i — Mo/2, the moment the planet is ejected. 
Similarly, using Eqs. (I39p and (|40p . one can show 
that aout = oo and adrc = ao(l-l-eo). Therefore, the 
circularisation semimajor axis is at most twice the initial 
semimajor axis. When a planet is circularised, it is done 
so only momentarily; it can only retain such an orbit if 
the mass loss is suddenly stopped at that moment. For 
any planet that has been circularised, one can show that 
the semimajor axis will subsequently evolve as: 

apost — circular — arunaway | / — 180° • V / 

Therefore, the semimajor axis evolves through the e = 
transition smoothly, without changing its functional 
form. 

We test the goodness of these analytical approxima- 
tions by considering a close-in planet (at ao = 2 AU) in 
the robustly runaway regime of a supernova. Consider a 
IOMq progenitor which expels a — 0.5MQ/hour of mass 
past the orbit the planet. Thus, * ^ 62.4. When /o = 0° 
and /o = 180°, Eqs. ((33, ((Ml) and (gT]) replicate the ec- 
centricity evolution. Therefore, we set /o — 20° in Fig. 
1101 to show the extent of the deviation from the analytic 
approximation. In the figure, the thin black dashed lines 
represent the analytic approximation, which mimics the 
true eccentricity evolution to within 10% for all values of 
eo. As predicted by Eq. (|39|l . all planets are ejected be- 
fore the star loses half of its mass (at 10 hrs) . In Fig. 1111 
we set /o = 178°, just 2 degrees off an exact match, be- 
cause such a deviation from the analytics is more drastic 
than for deviations of /o = 0° . As /o deviates from 180° , 
the eccentricity turns up sooner, and becomes less cir- 
cularised. The approximation will mimic the semimajor 
axis evolution until the point at which the planet would 
have been circularised had /o = 180° . Note that these cir- 
cularisation instances occur when a is less than twice its 
initial value, in conformity with a^rc = ao (1 -f- eo). The 
dots in the left panel indicate when this circularisation 
would have taken place, and show that the semimajor 
axis evolution is unaffected. As predicted by Eq. (|42p . no 
planets are ejected until at least half of the star's mass is 
lost, and the highest eccentricity planets are ejected only 
in the limit of the star losing all of its mass (at 20 hrs). 
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Figure 10. Analytic approximations (thick black nearly-solid foreground lines, from Eas. l37l and l43p to the e and a evolution from 
the numerical simulations (background-coloured dashed curves) in a robustly runaway regime for /o = 20°. The planet at ap = 2 
AU is experiencing supernova-like mass loss of « = O.5M0/hour from a, fig = IOMq star {^q Ri 62.4). The lines with increasing 
dash length represent cq values of 0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, and 0.9, respectively. The analytical approximation is best 
for eo = 0.01, and reproduces all the eo curves from the full numerical integrations to within 10%. All the planets are ejected 
before half of the star's mass is lost (see Eq. l39l l. 
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Figure 11. Analytic approximations (thin black nearly-solid foreground lines, from Eqs. 1381 1411 and 144 l l to the e and a evolution 
from the numerical simulations (background-coloured dashed curves) in the same runaway regime as in Fig. 1101 but for /o = 178°. 
In the limiting case of /q = 180°, the eccentricity decreases until reaching zero. The dots in the left panel indicate when this 
would have occurred (at ao [1 + eo]). None of the planets are ejected until at least half of the star's mass is lost, and the highest 
eccentricity planets are not ejected until almost all of the star's mass is lost (see Eq. l42l l. 



2.7 Impulsive Regime Evolution 

One may treat the entirety of stellar mass loss under 
the impulse approximation, when the mass is lost instan- 
taneously. This situation corresponds to 'i'o — > oo, an 
asymptotic runaway regime. Let the subscripts "i" and 
"f" represent the initial and final values, £ the (uncon- 
served) specific energy of the system, and r and v the 
position and velocity of the planet. Then 

Gui _ Gfj,i 



£. = ^v. 



(46) 



and 



1 2 



Gm/ 



>0 



2 ^ rf 

assuming that the planet is ejected. 



(47) 



Now assume /i/ = P^i, where < /? ^ 1. We can 
obtain a condition for ejection by eliminating Vi = Vf 
from the equations and setting n — rf. Doing so gives: 



P> 



1 + ep -I- 2eo cos /p 
2(l-|-eocos/p) 



(48) 



We illustrate the phase space suggested by Eq. (|48|l in 
Fig. 1121 Below each curve of a given ep, the planet is 
ejected. Note how the region around /o = 180° highlights 
a stable region, one for which the highest eccentricity 
planets are the most protected. This situation is reflected 
in the finite '^ runaway regime, and demonstrated in Figs. 
[9] and 1111 Although the highest eccentricity planets are 
the most protected at /o « 180° (apocenter), they are 
the least protected at /o ~ 0° (pericenter). The tendency 
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for nearly circular planets to be ejected is independent of 
true anomaly. 

One curious connection between the impulse approx- 
imation and the bifurcation point is that the two inflec- 
tion points of Eq. (|48|) satisfy Eq. ([23} . The value of /3 at 
these points, /?infl, is given by: 



Impulse Approximation 



Anfl ^ 



l + Seg 



(49) 



Therefore, 1/4 ^ /3infl Sj 1/2. These points are marked as 
dots in Fig. 1121 and are connected by the analytic curve 
from Eqs. ^ and (O. 

Returning to Eq. (|48p . note that in the limit of 
/ -> 0° and / -s> 180°, one recovers Eqs. ^ and ((42)) . 
Additionally, if the inequality in Eq. (fl5)) is solved for 
cos / and then bounded by its maximum value, then ejec- 
tion is impossible if 



/3> 



eo 



(50) 



This condition demonstrates that for an initially circular 
planet, at least half of the star's mass must be lost in or- 
der for there to be a possibility of ejection. For a highly 
eccentric planet, however, just a slight mass loss might 
be enough to eject it. Whether or not a planet's high ec- 
centricity serves as a protection mechanism is dependent 
on its / value, which relates how close the planet is to 
pericenter or apocenter. 

As an example, consider a circular ring of massless 
particles uniformly distributed in true anomaly at any 
separation from a star of any nonzero mass. If over half 
of the star's mass is lost instantaneously, then the entire 
particle ring will be ejected. Otherwise, all the particles 
will remain bound. Now consider an eccentric ring where 
all particles have e = 0.9. If the parent star instanta- 
neously loses 60% of its mass, then only « 11% of the 
ring will remain bound to the star. 



3 EXCITATION AND EJECTION IN 
REALISTIC SYSTEMS 

3.1 Overview 

We can now apply the theory developed in Section 2 to 
realistic systems. The field of stellar evolution is exten- 
sive and touches on several areas of astrophysics. We can- 
not hope to cover the entire phase space in detail in one 
paper. However, by focusing on a single phase of stel- 
lar evolution and considering constant mass loss in most 
cases, we will attempt to provide preliminary statistics 
and order-of-magnitude analysis for the entire progenitor 
stellar mass range up to ISOM©. We perform detailed 
nonlinear simulations only for the 2Mq ^ fJ^o ^ 7i\f0 
regime, whose stellar evolutionary tracks lend themselves 
well to this study. 

The evolution of stars is almost entirely determined 
by its ZAMS (Z ero- Age Main Seque nce) metallicity con- 
tent and mass (jWooslev et al.ll2002l \ These two factors 
determine how mass is lost later in life through winds. 
Because this correlation is so poorly known, mass loss 
prescriptions are often treated as a third independent 
parameter for tracing stellar evolution. To avoid detailed 
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Figure 12. Ejection prospects in the impulse approxima- 
tion ('I'o -^ oo). Plotted is the fraction of stellar mass re- 
tained (= P) vs. /o for 10 curves of increasing dash length for 
eo = (0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9). The planet 
remains bound in the regions above the curves and is ejected 
in the regions below the curves. The dots refer to the inflec- 
tion points of the curves, given by Eq. II49I I and which satisfy 
Eq. II23I I. Note how highly eccentric planets are especially pro- 
tected from ejection near apocenter (/ fa 180°), but arc prone 
to ejection near pericenter (/ fa 0°). 



modelling involving integration of the stellar evolution 
differential e quations, w e rely heavily on the empirical al- 
gebraic fits of lHurlev et a l. (2000) to model the evolution- 
ary tracks of stars of most mass, metallicity and mass loss 
rate properties. These evolutionary tracks demonstrate 
that mass loss i) can occur in multiple stellar phases, 
ii) is often prominent in just one stellar phase, and iii) 
is always monotonic but typically nonlinear in any given 
phase. All stellar evolutionary phase names used here are 
defined in their seminal work. 

We use t he m ass loss prescriptions provided in 
iHurlev et al.l l|2000f ). which inclu de the Reimers law 
on th e Red Giant Branch (RGB; iKudritzki fc ReimersI 
Il978l ). a stea dy superwind asym pt otic giant branch 
prescription l|Vassiliadis fc WoodI ligGS]). a high- 
mass loss prescription IjNieuwenhuiizen fc de Jageij 
199Q|), a Wolf-Ray et-like mass loss prescription 
(iHamann et al.l 119951) and a luminous blue variable 



law ( Humphreys fc DavidsonI [l99J) . The Reimers pre- 
scription is in particular widely used for giant branch 
evolution, and beyond: 



dt 



= jy (4 X 10" 



L-,(t)n-,(t) Mq 



M* 



yr 



(51) 



where L* and 7?.* are the stellar luminosity and ra- 
dius and r; is a dimensionless co efficient. We adopt 
the "typical" value fo r r) of 0.5 (|Hurlev et all l2000l : 
ISchroder fc Cuntzll2005l '). 

We divide the stellar mass phase space into 5 
regimes, which are approximately separated at IM©, 
2Mq, 7Mq, and 2OM0, based on stellar evolutionary 
properties. 
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3.2 Numerics and Checks 

Although certain regimes of evolution can be modelled 
well by algebraic formulas, the lack of a complete closed- 
form analytical solution to Eqs. ((3])-(l7|) suggests that nu- 
merical integrations are needed to model the evolution 
of the full two-body problem with mass loss. We evolve 
planetary orbits in this section by using numerical inte- 
grations. 

In these integrations, one may incorporate mass loss 
i) as a separate differential equation, ii) by explicitly re- 
moving mass from the primary according to a given pre- 
scriptio n, or, alternatively iii) by ad ding mass to the sec- 
ondary (JDebes fc SigurdssonI I2OO2I ) . As a check on our 
results, we have reproduced each curve in the breaking 
of adiabaticity regime in the left panel of Fig. [3] ('I'o ~ 
0.023) with both i) integration of the orbital elements 
(Eqs. I3I7|I plus integration of a separate mass loss differ- 
ential equation in the Mathematica software program, for 
13 digits of accuracy and precision and with a working 
precision equal to machine precision, and ii) integration 
of the Cartesian equations of motion w ith the hybrid in - 
tegrator of the N-body code. Mercury (|Chambersll 19991 ). 
with a maximum timestep of 1 yr and with mass explic- 
itly being removed from the primary at each timestep. 

However, we warn future investigators that in sys- 
tems which ultimately do not obey the adiabatic ap- 
proximation, the dynamical evolution is sensitive to the 
evolution of /. Therefore, in numerical integrations, par- 
ticularly for nonlinear mass loss prescriptions, how one 
discretises the continuous mass-loss process can quali- 
tatively affect the resulting evolution. A discretation of 
mass loss is mimicked in reality by instantaneous bursts 
of primary mass lost beyond the orbit of the secondary. 
Therefore, a detailed study of an individual system with 
a given mass loss prescription will require a numerical 
integration where the time between discrete decreases in 
the primary mass should be less than the (largely un- 
known) timespan of discontinuous patterns in the mass 
loss modeled. Here, we just seek to demonstrate the in- 
stability in the general two-body mass loss problem and 
achieve representative statistics on ensembles of systems. 
In the 2Mq ^ a^o ^ 7Mq regime, which features nonlin- 
ear mass loss, we model planets with 50 AU ^ ao ^ 10"" 
AU. Therefore, we set a maximum possible timespan of 
1 yr (the same value used to reproduce Fig. ^ between 
mass lost; in order to achieve mass loss on this scale, 
we interpolate linearly between the outputs from the 
largely nonlinear st ellar evolutionary track outputs from 
iHurlev et al.l (l200d). We th en run the simulations with 
Mercury's ( Chambersll 19991 ) hybrid integrator. 



3.3 The Stellar Mass Spectrum 

3.3.1 The fio < IMq regime 

Sub-solar mass stars experience quiescent deaths, some 
of which are theorized to last longer than the age of 
the universe. However, stellar tracks computed from the 
iHurlev et al.l (2000) code indicate that the most massive 
members of this group (/lo > 0.7Mq) may pass through 
multiple stages of evolution, and eject up to half of their 



initial mass in the Red Giant Branch (RGB) stage. Low 
metallicity /lo ~ O.SMq and /lo ~ O.9M0 stars do so 
on the RGB over ~ 100 - 200 Myr. If this mass is lost 
uniformly, then ^0 ~ 0.011, meaning that the system is 
likely to start losing its adiabatic properties. Simulations 
of constant mass loss confirm that the change of the ec- 
centricity of an Oort Cloud at ao = 10^ AU will vary 
from ~ 0.01 (for particles with eo — 0.90) to ~ 0.1 (for 
particles with eo = 0.01). The mass loss is not strong 
and quick enough to eject the particles, and objects with 
semimajor axis less than ~ 10* AU (which would yield 
* $; 0.00036) are robustly in the adiabatic regime. This 
regime of the motion might change, however, due to non- 
linear modeling. This might reveal short bursts of mass 
loss causing '^ to increase sharply over the corresponding 
burst timescale. 



3.3.2 The IMq < /lo < 2Mq regime 

Roughly half of all known planet-hosting stars, including 
the Sun, lie in this progenitor mass regime, motivating 
detailed analyses of these systems. We defer such analy- 
ses to future studies because of the complex multi-phasic 
evolutionary path these stars are prone to follow. 

As an example, assuming rj = 0.5, the Sun will even- 
tually lose a total of 48% of its original mass: 24% dur- 
ing the RGB, 4% during core-He burning, 13% during 
the Early Asymptotic Giant Branch (EAGB), and 7% 
during the Thermally Pulsing Asymptotic Giant Branch 
(TPAGB). All these phases of mass loss are nonlinear and 
occur on different timescales. If instead rj = 0.3, then the 
mass loss percentages will change drastically: 13% dur- 
ing the RGB, 2% during core-He burning, 4% during the 
EAGB, and 28% during the TPAGB. Other examples 
show that slightly increasing the progenitor mass from 
I.IMq to 1.2Mq can have a similarly large effect on what 
mass is lost when. 

We can, however, provide some rough estimates 
of planetary evolution through representative numerical 
simulations assuming constant mass loss over one phase. 
Stars in this mass regime may lose over 60% of their orig- 
inal mass, most of which either in the RGB (particularly 
for values of r; ^ 0.8) or the TPAGB (for lower rj and 
po > I.3M0). The duration of RGB phases for these 
masses are ~ 100 Myr, and will yield only minor eccen- 
tricity increases at a = 10^ similar to those from sub- 
Solar masses. However, the duration of TPAGB phases 
in this mass regime is ~ 0.1 — 1.0 Myr. Constant mass 
loss over this period of time for ^0 = 1.0A/q — 1.3Af0 
can cause up to 20% of an Oort Cloud at 10^ AU 
(^ Ri 3.0) to be ejected, and raise the eccentricity of 
an initially circular planet at 10* AU {'i' ~ 0.096) to 
~ 0.25. We obtained these figures by sampling 8 evenly 
spaced values of /o for each of the following 10 values 
of eo: 0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9. This ef- 
fect is pronounced with progenitor masses approaching 
2M0 and losing up to 70% of their initial mass. 

Therefore, Oort clouds are in jeopardy of partially 
escaping or being moderately disrupted in systems with 
similar progenitor masses to the Sun. The comets cannot, 
however, drift into the inner regions of the system (see 
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Ea. l21|) . The widest-orbit planets at ~ 10* AU may expe- 
rience a moderate eccentricity change of a few tenths, and 
might be ejected depending on the nonlinear character of 
the mass loss. Future multi-phasic nonlinear modelling 
will better quantify and constrain these effects. 



3.3.3 The 2Mq sJ mo < TM© regime 

This mass regime is well-suited for this study because 
here, ~ 70%— 100% of a star's mass loss occurs in a single 
phase, the TPAGB, regardless of the values of r), [Fe/H], 
or ^Q. Therefore, by modelling the nonlinear mass loss in 
this one phase, we can make definitive conclusions about 
this region of phase space. Additionally, the duration of 
this phase is short, typically under 2 Myr, and therefore 
feasible for numerical integration of planets at distances 
of just a few tens of AU. 

We consider two progenitor star metallicities, a 
"low" metallicity ([Fe/H] = 0.0001), and Solar metal- 
licity ([Fe/H] = [Fe/H]© = 0.02), both with ry = 0.5. In 
the low metallicity case, we utilize 9 TPAGB evolutionary 
tracks that range from /io = 2Mq—6Mq, in increments of 
O.5M0. In the solar metallicity case, we utilize 13 TPAGB 
evolutionary tracks that range from no — 2Mq — 8Mq , in 
increments of O.SA/©. Beyond these upper mass limits, a 
star would undergo supernova for the stated metallicities. 
The evolutionary tracks are plotted in Fig. 1131 Note that 
the initial masses indicated on the plots do not exactly 
represent no; the small (< 10%) mass loss which occurred 
between the main sequence and the start of the TPAGB, 
typically in the Core-He burning and EAGB phases, has 
already been subtracted. For all of the tracks except the 
low metallicity /xo = 2Mq track, most of the mass loss 
occurs within a short 10* yr scale indicated by the sharp 
downturn in the curves. However, note that between the 
start of the TPAGB phase to this intense mass loss pe- 
riod, over a period of ~ 0.7— 1.5 Myr, the stars typically 
lose ~ O.5M0 worth of mass. After the intense mass loss 
burst, effectively no more mass is lost from the system. 
Integrations for the two lowest-mass tracks for [Fe/H] — 
[Fe/H]Q = 0.02 were begun 5 x 10^ yr after the start of 
the TPAGB in order to include the sharp mass loss fea- 
ture and consistently integrate all systems over the same 
period of time. 

For each of the 22 evolutionary tracks, we mod- 
eled 1200 planets as test particles and integrated the 
systems for 1.6 Myr, longer than the duration of the 
TPAGB phase for nearly all of the stellar tracks. The 
planets were all given randomly chosen values of the 
initial mean anomaly, and were split into 8 groups of 
150. Each group of planets was assigned an ao value of 
50, 100, 500, 1 X 10^ 5 X 10^ 1 X 10*, 5 x 10*, and 1 x 10^ 
AU. Each group of 150 planets was split into three sub- 
groups of 50, each of which was assigned an eo value of 
0.01, 0.5 and 0.9. 

We compute the percentage of each group of 150 sim- 
ulations of a given semimajor axis and initial progenitor 
mass which become unstable. We define instability by 
whether or not the planetary eccentricity reaches unity. 
Figure [T?] reports the results. Because of the nonlinear 
nature of the mass loss, here our mass loss index from 



Eq. (|15l) breaks down. However, we can say roughly that 
the duration of the greatest mass loss is comparable to 
a planet's period at 500 AU (the blue curves with dia- 
monds). At approximately this semimajor axis we expect 
a planet to be in the transition region between adiabatic- 
ity and runaway. This curves on the plot qualitatively 
corroborate this expectation: orbits tighter than 500 AU 
are stable and adiabatic, orbits wider than 500 AU are 
unstable and runaway, and orbits at 500 AU are a bit 
of both. Most of the planets in the widest orbits become 
unstable, but they cannot all become unstable for a large 
enough sample of randomly chosen values of /o because 
some of these values will be close to 180°. As demon- 
strated by Fig. 111! in the high (e.g., 0.9) eo case, planets 
with /() close to 180° will be ejected only if the parent 
star loses over ~ 95% of its mass, a largely unrealistic 
scenario for any progenitor mass. Further, for the sim- 
ulations in Fig. 1141 note that beyond 1000 AU - in the 
robustly runaway regime - there is little correlation with 
instability percentage and ao. Equations (|37p . H38p and 
(|41|) help show why: at least for values of /o close to 0° 
and 180°, the eccentricity evolution is independent of ao. 
For the planets which remain bound, we consider the 
extent of their eccentricity excitation. Figure[TS]plots the 
eccentricity range experienced by bound planets averaged 
over all simulations with the same values of /io, ao, and eo 
but with different values of /o. The panels show that the 
eccentricity of the remaining bound planets for ao ^ 500 
AU is significantly excited (by several tenths). The ec- 
centricity of planets at ao = 50 and 100 AU on average 
can vary by a few hundredths, and 0.1, respectively. If a 
symbol in the legend does not appear on the correspond- 
ing plot, then no planets at that semimajor axis remained 
bound. The top two panels (eo = 0.01) exhibit a dearth 
of these symbols, a result one might expect from Fig. 1121 
If that figure is qualitatively representative of the situa- 
tion here, amidst strong nonlinear mass loss, then there is 
no value of /o which affords the lowest eo planets protec- 
tion. The horizontal lines on the middle two and bottom 
two panels of Fig. [15] display the value of 1 — eo; symbols 
above these lines indicate that the corresponding systems 
on average experience a net eccentricity decrease. These 
systems are more likely to be left with a planet whose or- 
bit is less eccentric than eo when mass loss is terminated. 



3.3.4 The 7Mq < /lo < 20Mq regime 

Generally, Solar-metallicity stars with 8Mq ^ no ^ 
20Mq are thought to undergo supernova and produce 
a neutron star. However, these bounds are approximate. 
Additionally, lower metallicity stars can begin neutron 
star formation and black hole formation at different 
values; repre sentative ones might be GM© and 18Mp, , 
respectively (JHeger et al.l |2003| ; lEldridge fc ToutI |2004| : 
iBelczvnski et al.ll2O10l ). These stars may eject ~ 50% — 
95% of their initial mass, most of which is in the su- 
pernova (ISmartt et al.ll2009l ) . Additionally, the minimum 
and maximum possible masses of the remnant neutron 
stars are constrained by physical principles. Typically ac- 
cepted values for the min imum and ma ximum are ~ IMq 
(jStrobel fc Weigell I2001I ) and 3Mq (|Kalogera fc Bavml 
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TPAGB Evolution for [Fe/H] = 0.0001 
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Figure 13. Thermally Pulsing Asymptotic Giant Branch (TPAGB) evolution for stars of low- and Solar-metallicities. Each colour 
represents a different evolutionary track. Initial TPAGB mass can be read off from the Y-axis. The four highest-mass gray curves 
for Solar metallicities were not computed for the low metallicity case because those stars would have undergone supernova. 
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Figure 14. Planetary ejection prospects for massive stars from 2Mq — SMq . Each data point is averaged over the 150 randomly 
chosen mean anomaly values and 3 selected eo values for each ao. The black fiUed-circle curves for 50 AU are hidden behind the 
100 AU curves. The ag = 50 and ag = 100 AU systems, which are in the adiabatic regime, remain bound. The oq ^ 10^ AU 
systems, in the runaway regime, become largely unstable. 



ll99d ): IClark et al.l ([2003) presents observational evidence 
for the upper bound. Therefore, this mass regime of stel- 
lar evolution is relatively well-constrained, and due to the 
nearly instantaneous mass loss, is very well suited for this 
study. 

The sudden nature of the supernova, combined with 
the great extent of mass lost compared to /io, place 
any orbiting planet immediately in the runaway regime. 
Therefore, we seek to determine what planets, if any, can 
survive a supernova. We hence choose parameters that 
favor survival, to see if this situation is possible. First, 
we s elect the minimum p ossible ao. Evolutionary tracks 
from lHurlev et al.l (l2000|) indicate that the minimum ex- 
tent of the pre-supernova stellar envelope (including pre- 



supernova mass loss) is about 
ao = 2 AU. 



2 AU, so we choose 



If the mass ejected from a supernova is considered 
to be isotropic, then this mass will collide with any or- 
biting bodies. This collision is likely to destroy smaller 
bodies. Large and/or massive planets, however, may sur- 
vive. Those planets which do survive might accrete some 
of the mass from the ejecta. Although doing so will cause 
a to decrease, this contribution, even at 2 AU, is negligi- 
ble compared to the a increase from all the other ejecta 
that is being blown past the planet's orbit. We are con- 
cerned with the amount of time the mass takes to pass the 
diameter of the planet. We can model mass loss in these 
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Figure 15. Eccentricity excitation of planets which remain bound during massive star evolution, for 2Mq — 8Mq. The left and 
right panels are for low metallicity and Solar metallicity, respectively. The top, middle and bottom panels are for eo = 0.01,0.5, 
and e = 0.9, respectively. Each data point is averaged over the 50 values of the mean anomaly sampled for the given fiQ, ag, and 
eo values. If no symbol is displayed, then none of the corresponding systems were stable. The horizontal lines indicate values of 
1 — eo; symbols above this line experience a net eccentricity decrease. 



systems by assuming an ejection velocity and a planetary 
diameter. 

Observations help constrain the velocity of this 
ejecta. Some diverse examples fo r differ ent types of Su- 
pernovae include: i) iFesen et al.l (2007) report Hubble 
Space Telescope observations which indicate that the 120 
yr average expansio n velocity of SN1885 is 1.24 x lO** ± 
1.4 X 10^ km/s, ii) iMazzali et"aLl l|2010l ) model spectra 



of SN2007gr, and find that the inner IMq of material 
is being ejec t ed at a velocity of 4.5 x 10'^ km/s, and iii) 
ISzalai et al.l (|201ll ) find that the maximum velocity of 
supernova ejecta of 2004dj during the nebular phase is 
approximately 3.25 x 10^ km/s. One theoretical investi- 
gation claims that ejecta velocity c an reach 2 x lO** km/s - 
3 X 10* km/s IjWooslev et al.ll 19931 ) , and another demon- 
strates that (surface) piston speeds of 1 x 10* km/s - 
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2 X 10'' km/s "covers the extremes from a sudden (en- 
ergy deposition over Is) to a slow- developing explosio n 
(energy deposition over ~ 100 ms)" IjDessart et al.ll2010f ). 
As exemplified by these examples as w ell as the compila- 
tion in Fig. 1 of lHamuv fc Pintd (|2002r ). a typical range is 

V = 10^ — 10* km/s; let us then assume the lower bound 

V = 10^ km/s. 

Further, the highest known exoplanet radius is less 
than twice Jupiter's radiu^jO, so let us assume this value 
for our planet. We can then test the extremes of the total 
mass lost (= Meje) based on the progenitor mass and 
remnant mass bounds. 

We assume the mass is blown past the or- 
bit of the planet isotropically and at a constant 
value, and we consider 36 uniformly distributed val- 
ues of /o. For each, we adopt 10 values of bq 
(0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9). We simulate 
these 360 systems in each of four scenarios: with the two 
extreme values of /io {GMq and 20Mq) and two extreme 
values of the remnant mass {IMq and 3Mq). For all these 
cases, ^^o ~ 1 — 2, placing these systems in the weak run- 
away regime at i = 0. The reason why ^I'o is not higher 
for such a great mass loss rate is because ao is so low 
(2 AU). However, to determine the endpoint of orbital 
evolution, one needs to combine an estimate of ^o with a 
mass loss duration time (or a remnant mass, for constant 
mass loss), which is independent of ^o- This is why the 
endstates can change drastically for two systems even if 
their initial mass loss indicies are equivalent. 

Figure [16] displays the result of our simulations. The 
figure demonstrates that an appreciable number of plan- 
ets can survive, but only in the extreme case of the su- 
pernova ejecting just half of the progenitor mass, and 
only if /o is closer to 180° than to 0° . In the more realis- 
tic cases of greater mass loss during supernova, the only 
planets which may survive must have /o « 180°. This 
initial condition appears to be the only protection mech- 
anism against ejection for robustly runaway (see Fig. lllf) 
or impulsive (see Fig. I12p systems which lose most of 
their mass. The impulsive limit can further help explain 
Fig. 1 161 through Eq. (|48|l : the top, black curve with open 
circles corresponds to /3 = 1/2, and the highest initial 
eccentricity we sampled in the simulations was eo — 0.9. 
Therefore, Eq. (|48|) gives cos /o < —0.9, meaning that all 
planets with 154° < /o ^ 206° should remain stable. The 
numerical simulations confirm the theory. Additionally, 
Fig. [12] confirms why no planets with initial true anoma- 
lies within 80° of pericenter survive, even when just half 
of the star's mass is lost. 

In fact, if we decrease or increase the mass loss rate 
(and hence 'I'o) by an order of magnitude (to either 

V = 100 km/s or u = lO* km/s), and rerun our simu- 
lations, we reproduce Fig. [16] closely. Therefore, in such 
runaway regimes, the evolution becomes independent of 
the mass loss rate above a certain critical mass loss rate. 
The realistic implication of this finding is that the par- 
ticular choice of ejecta velocity assumed for a supernova 



Surviving Supernovae at ao = 2 AU 
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Figure 16. Survivability of a tight-orbit (2 AU) planet during 
supernova. The progenitor mass and ejected mass is given by 
Mq and Mcjc, respectively. Each data point is based on 10 sys- 
tems with eo = (0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9). For 
all cases, \&o ~ 1 — 2. The plot demonstrates that the only way 
planets may remain bound after a supernova blast is by ini- 
tially residing in a narrow region of true anomaly space. This 
behavior was predicted in Sections 12.6.21 and 12.71 and specif- 
ically in Fig. 1121 Equation II48I I demonstrates why all planets 
on the black curve with open circles with 160° ^ /o ^ 200° 
survive, and why all planets near pericenter are ejected. Al- 
most all planets that withstand a supernova are ejected unless 
the percentage of the progenitor's mass lost is the theoreti- 
cal lower bound (sd 50%) for supernovae of fa 7Mq — 20Mq 
progenitor masses. 



is unimportant, as long it is assumed to be higher than a 
critical minimum value. 



3.3.5 The fio > 20Mq regime 

There is great uncertainty regarding how the highest- 
mass stars lose mass and in what amounts. The possi- 
bilities for planetary evolution around these stars are in- 
triguing, and can be simulated once a model has been 
adopted for a particular star. Stars in this regi me gen- 
erally beco me neutron stars o r blac k hole^j. Hege r et al.l 
l|2003l ) and [Eldridge fc Tout! (|2004l ) outline these poten- 
tial stellar fates as a function of initial progenitor mass 
and metallicity. 

Black holes may form with or without a supernova. 
In the latter case, mass is still lost during core collapse. 
Quantifying the extent and timescale of this mass loss is 
crucial for determining the fate of any orbiting planets. 
This process is thoug ht to last on the orde r of tenths 
of seconds to seconds (JQ'Connor fc Ottll201ll ). The mass 



^ http://exoplanet.eu/ 
■* http://exoplanets.org/ 



^ In rare cases, at the very lowest metallicities, pair instability 
supernovae will destroy the entire star and leave no stellar 
remnant. 
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lost during this process has been m odeled to be as much 
as 1-2 Mq (IBe czvnsk et al.ll2010l l. However, this value 
could be zero; iFrveiJ ( 1993 ) argues that for progenitor 
masses above AQMq, the fi nal black ho l e mas s could be 
as large as the progenitor. IZhang et al.1 ([2003) indicates 
that stars at even greater masses, with lOOM© < /xq < 
26OM0, may explode completely and leave no remnant. 
These theoretical treatments are poorly constrained by 
observations. Howe ver, observations do su ggest that stars 
up to 3OOM0 exist (jCrowther et al.ll2010l ). 



Hence, unlike in the previous subsections, stellar 
mass evolution here remains qualitatively uncertain, as 
the amount of initial mass lost could be any value up 
to 100%. Therefore, we can frame our cursory analysis 
in this section by considering what percent of the pro- 
genitor's mass must he lost in order to produce ejec- 
tion or excitation. The large progenitor masses in this 
regime promote adiabaticity, as indicated by Eq. (|15p . 
and hamper prospects for planetary ejection, as indicated 
by Eq. (|36p . However, mass could be lost through super- 
winds at a g reat rat e of 10~*M©/yr l|Dessart et al.ll2010l : 
lYoon fc Ca nticUo 2 0101 ). which might offset the stabilis- 
ing effect of the large magnitude of the progenitor mass. 



We can provide a preliminary overview of the impact 
large progenitor masses with superwinds would have on 
the survivability of planets. We consider three (strong) 
mass loss rates, a = {10~*, 10~^, lO~®}M0/yr and pro- 
genitor masses up to ISOMq. For all these cases, unless 
'^o <; 10* AU, the planetary evolution will be primarily 
adiabatic, as ^^o <C 1. We find that at least ~ 80% of a 
progenitor's mass must be lost for any planet at ao ~ 10'^ 
AU to be ejected by any of these winds. However, for 
planets at ao ~ 10^ AU, a mass loss of a = lO~''M0/yr 
does place the planet in a runaway regime. In this regime, 
for /o = 0°, the progenitor needs to lose just a few per- 
cent of its initial mass to eject the highest eccentricity 
planet, and roughly 50% of its mass to eject initially 
circular planets. These results conform to expectation 
from Eq. (|39p . and hold for all progenitor masses from 
2OM0 — 15OM0. Therefore, even without appealing to 
core collapse or weak supernova, the mass loss from the 
highest mass stars in the universe can blow away any 
remaining Oort Clouds. 



Detailed modelling of secondaries evolving amidst 
the complex evolution of stars in this mass regime is a 
ripe topic for future studies. Although the mass lost in 
core collapse can approach zero, the nearly-instantaneous 
timescale for the mass loss might have a sudden pro- 
nounced effect on the planetary orbit. Further, fallback 
of mass from a weak supernova explos ion onto a neutron 
star lasting "seconds to tens of hours" IjHeger et al.ll2003r ) 
could trigger a black hole. This fallback will cause a still- 
bound planet's semimajor axis to decrease. Also, for stars 
that explode away almost 100% of their mass, one may 
investigate the minimum amount of mass that could re- 
main and still bind a planet. In this case, the planet's 
mass will become important. 



4 DISCUSSION 



4.1 Oort Clouds 



No Oort Clouds have been observed. However, comets 
thought to originate in the Sun's Oort cloud have been 
observed, and have motivated several stu dies which es- 
timate the orbital extent of these bodies. iLevison et al.l 
(.20101) claims the Oort cloud extends to ~ 10^ AU, and 
iDvbczviisk i (2002) claims that this is a "typical" value for 
the outer boundary. Although planetary material might 
exist throughout the sca ttered disk from the Kuiper Belt 
to the Oort cloud (e.g. iLeto et al.l 



(JDuncan et al.lll987l : [Gardner et al.l 



200s ) , some authors 



2OIII ) have set an in- 



ner boundary at several thousand AU. Other studies fo- 
cus on a supposed break in the Oort cloud, separating 
it into an "inner" and "outer" region. This bifurcation 
is claimed to occur at at ~ 2 x 10"^ AU ()Hillsl Il98ll : 
iKaib fc Ouinnll2008l : iBrasser et aIll201G[ ). 

These estimates pertain to the Solar System only. 
Oort clouds around other stars may exist. Stars born 
in more dense clusters will have more comets deposited 
into their clouds than did the Sun. lKaib fc QuinnI (|2008h 
simulate four different primordial environments (with no 
cluster, and three clusters with densities of 10, 30, and 
100 stars per cubic parsec) and find that all produce sim- 
ilar "outer" (a > 2 X lO'' AU) Oor t clouds and qualita - 
tively different inner ones. Further. iBrasser et al.l (|201Gr ) 
consider the different types of Oort clouds which may be 
formed around other stars as a function of galactocentric 
distance. At large galactocentric distances (> 14 kpc), 
they find that some (> 10%) Oort cloud constituents or- 
bit beyond lO'^ AU. 

All these estimates suggest that the majority of stel- 
lar mass progenitors, including the Sun and those of sub- 
solar mass, will excite the eccentricity of Oort Clouds dur- 
ing stellar evolution. Most of these Oort Clouds will lose 
material to interstellar space. Assuming that the comets 
are roughly distributed uniformly in true anomaly, then 
only a fraction will survive. This fraction is highly de- 
pendent on the duration of mass loss. The remaining 
comets will as s ume a differential eccentricity distribution. 
IBrasser et al.l (J2010l l focus on galactic tides and how they 
strip off Oort cloud constituents. Indeed, Oort clouds may 
not even survive to the post main-sequence phase. If they 
do, the Galactic tide will be stronger relative to the star's 
gravity for any surviving Oort cloud objects, so the stable 
region that the Oort cloud can occupy will have shrunk 
at the same time that the bodies' orbits are expanding, 
potentially leading to even more ejections. Further, as a 
star loses mass, its gravitational influence within its stel- 
lar neighborhood will shrink and be encroached by the 
potential wells of stellar neighbors. 

However, as demonstrated by Fig. 3 of lHiguchi et al.l 

(120071 ) ■ galactic tides often need Gyr of evolution in or- 
der to cause an appreciable change of a comet's orbital 
elements. Short-lived massive stars won't often provide 
galactic tides with this opportunity before stellar mass 
loss becomes the dominant perturbation on the comets. 

More detailed modeling of Oort clouds could enable 
investigators to link mass loss from a white dwarf progen- 
itor with the cometary population of the resulting white 
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dwar f (see lAlcock at al.1 llQSd and iParriott fc Alcockl 
119981 ). Additionally, one should also consider the differ- 
ence in the stellar wind velocity at Oort Cloud distances 
versus its escape velocity when it lea ves the star. As ob- 
served bv lDebes fc SigurdssonI (J2002l ). because the wind 
crossing time is typically longer than the Oort Cloud or- 
bital timescale, winds which have slowed will enhance the 
system's adiabaticity. 



stability format ion models can fully form p lanets be- 
yond ~100 AU iDodson-Robinson et ai]|2009l '). Embryos 
and/or partially-formed planets that were scattered be- 
yond ~ 10"^ — 10^ AU will undergo the same dynamical 
evolution due to stellar mass loss as a fully- formed planet. 
This situation might arise around short-lived, high-mass 
stars, where the timescale for core accretion might be 
longer than the mass loss timescale. 



4.2 Wide-orbit Planets 

Initially, exoplanet discovery techniques were not well- 
suited for detecting planets which reside beyond ~ 6 
AU on decade-long timescales, and this region remained 
relatively unexplored until the mid-2000s. However, new 
observational techniques and carefully targeted surveys 
are increasing t he likelihood of uncove ring planets on 
wide orbits fe.g. ICrepp fc Johnsonll201ll ). The discover- 
ies of the fou r planets with a ~ 15, 24 , 38, 68 AU orbit- 
ing HR 8799 (JMarois et al.ll2008l . [201ol') and the a kj 115 
AU planet orbiting Formalhaut (JKalas et al.l |2008|) re- 
vealed that wide-orbit (a > 10 — 100 AU) planets do 
exist and incited great interest in their formation and 
evolution. Additionally, at least 10 wider-orbit compan- 
ions which may be massive planets that are close to the 
brown dwarf mass limit have been det ected. Like For- 
malh aut b, the companion to GQ Lup IjGuenther et al.l 
I2OO5I ) is thought to satisfy 100 AU < a < 20 AU. 
Cor npanions around A B Pic (jChauvin et al.ll2005h. Qph 
11 llClose et all |2007| ) and CHXR 73 (JLuhman et all 
120061 ) all harbor semimajor axis betwee n 200 AU and 
300 AU, and those orbiting C T Cha dSchmidt et"aLl 
I2OO8I ). IRXS J160929.1 -210524 (iLafreniere et al.l |201G| ) 
and CSC 06214-00210 l|lreland et al.1 I2OIII ) satisfv 300 
AU < a < 500 AU. Companions with 500 AU < a < 1000 
AU i nclude those orbiting UScoCTIO 10 8 l|Beiar et al.l 
I2OO8I). HIP 78530 jLafreniere et al.ll201ll ) and HN Peg 
B ( Leggett et al.ll2008l ). The three potentially planetary 
companio ns with the widest kno wn orbits are SR 12 C 
fllOO AU.lKuzuhara etHlbOllh . Ross 458 b (1168 AU, 
Goldman et alll2010l ). and WD 0806-661B b (2500 AU, 



Luhman et al.ll201lh . Theoretical models place the mass 



of the a — 2500 AU object at 7 Jupiter masses. Our study 
is particularly relevant to such wide-orbit companions. 

During stellar evolution, wide-orbit planets in isola- 
tion will behave equivalently to Oort cloud comets, and 
can be ejected with similar ease. Planets may be mutually 
scattered out to distances of ~ 1 AU while remaining 
bound to their parent systems (| Veras et al.l 120 09); be- 
yond this distance, over time the effects of passing stars 
are likely to strip the planet from the system. Another 
mechanism for producing wide-orbit planets is capture 
from other stars, or passing free-floaters. There is still a 
possibility that the Sun contains a massive, very- wide or- 
bit companion. iFernandea (|201ll ) discusses the prospects 
for detecting a wide-orb it (> 10'* AU) Jovian mass com- 
panion to the Sun, and lMatese fc Whitmird (|201ll ) sug- 
gest that there is evidence for such a companion re- 
siding in the Sun's outer Oort cloud. Regardless, such 
planets are very unlikely to have formed in these en- 
vironments; neither core accretion nor gravitational in- 



4.3 Multiple Planets 

Introducing additional bodies in the system, such as a 
second planet, or a belt of material, could si g nifica ntly 
complicate the evolution. iDebes fc SigurdssonI 120021) in- 
vestigate the first scenario, and JBonsor et al.l ((2011) the 
second. In both cases, the characteristics of their N-body 
simulations demonstrated that the systems they studied 
were in the adiabatic regime. In this regime, where stellar 
mass loss produces quiescent adiabatic eccentricity exci- 
tation on the order of $0 (see Eq. I17[) . the eccentricity 
variation of the second planet or belt particles can then 
be attributed solely to the other planet. Additionally, the 
orbit of the true anomaly is only negligibly affected by 
mass loss in the adiabatic limit. Thus, the main contri- 
bution of the stellar mass loss in their studies is through 
the well-defined (Eq. [2| increase in semimajor axis of all 
objects in the system. 

Including additional planets in situations where \&bif 
is reached and/or exceeded represents several of the nu- 
merous potential extensions to this work. The frequency 
of planet-planet scattering and the resulting free-fioating 
planet population in the midst of semimajor axis and ec- 
centricity variations from stellar mass loss are important 
issues to be addressed. Other situations to consider are 
how planets may stay locked into or be broken from sec- 
ular and mean motion resonances, and how instability 
timescales are affected. 



4.4 Ftee-Floating Planets 

The ejection of planetary material, whether it be in the 
form of partially-formed planets, fully-formed planets, or 
comets, might contribute to the free-fioating mass present 
and potentially detectable around dead stars. Evidence 
for t he existence of free-fio a ting planets has been mount- 
ing (JLucas fc Roche _ 200d : IZapatero Osorio et al.l I2OO0I , 
|2002| : JBihain et al.ll2009l ) and was recently highlighted by 
a re port of potential detections of 10 free-fioating plan- 
ets (jSumi et al.ll201ll ). Also, the capability may exist to 
distinguish between free-fioaters and bound wid e-orb it 
planets up to semimajor axes of « 100 AU (iHani: 2006l ). 

Assuming that the same amount of planetary mate- 
rial was distributed equally among stars of all progenitor 
masses, then « 7Mq — 2QMq progenitors are by far the 
most likely stars to produce free-fioating materialj, fol- 
lowed by stars in the ~ AMq — ?>Mq progenitor mass 



^ One potential indication of tlie origin of supernova- produced 
free-floaters is tlieir space velocities; neutron star "kicks" cause 
the true space velocities of yo ung pulsars to refle ct tlie (high) 
speed of the supernova ejecta IIHobbs et al.ll2005l) . 
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range (see Fig. I14|) . The ability of stars with fio > 20Mq 
to produce free floating material is unclear and is largely 
dependent on the evolutionary models used. For a given 
progenitor mass, metal-poor and/or metal-rich stars may 
be prone to ejecting planets. However, because metal- 
rich stars are slig htly more likely to ha rbor planets than 
metal-poor stars (jSetiawan et al.ll2010f ). the metal-poor 
stars which are dynamically prone to planetary excita- 
tion might not initially harbor planets. 

Detailed modelling of the galaxy's free-floating 
planet population requires 1) an initial mass function, 
2) better statistics for planets orbiting stars other than 
Sun-like hosts, 3) knowledge of how many planets inhabit 
wide orbits at for example, a = 10'^"^ AU, and 4) better 
knowledge of the ability for ~ IMq — 2Mq stellar-mass 
progenitors to eject planets. Depending on these results, 
stellar evolution might be the primary source of free- 
floating planets. Alternatively, if, for example, a negligi- 
ble number of planets are shown to inhabit orbits beyond 
a = 10'' AU, then the dominant source of free-floating 
planets would likely lie elsewhere. 



4.5 Pulsar Planets 

Our results suggest that very few first-generation pulsar 
planets exist. Such planets would have had to reside far 
enough away from the expanding progenitor envelope to 
not be disrupted pre-supernova, and then survive the su- 
pernova. Assuming a uniform distribution of true anoma- 
lies, only (180° -/crit)/180° « 11% of planets would have 
a fighting chance to survive due to the additional time 
they would take to initially decrease their eccentricities. 
Even then, their initial eccentricities would have to be 
high enough, and the mass loss duration short enough, to 
outlast the supernova. These results suggest that unless 
pulsars can readily form planets or capture them from 
other systems, pulsar planets should be relatively rare. 

However, if the pulsar planet survived engulfment 
from the expanding pre-supernova stellar envelope, then 
its semimajor axis might be small enough to remain 
bound during the supernov a. There is one planet, HIP 
13044 (jSetiawan et al.ll2010^ ■ who potentially could have 
survi ved residing inside its star's envelope (jBear et al.l 
I2OIII ). The spiral-in time of the planet could have ex- 
ceeded the short duration (~ 100 yr) of the RGB expan- 
sion and engulfment, allowing the planet to survive. If 
close-in (< 1 AU) pulsar planets survive in a similar way, 
their final eccentricities could be any value (see Fig. Ill|l 
but their semimajor axes will have increased by many 
factors. The three pulsar planets orbiting PSR1257-I-12 
(jWolszczan fc Fraill1l99^ : IWolszczanlll994r) all have a < 
0.5 AU. If they are first-generation planets, then ao < 0.1 
would have held true for each. At such a small semimajor 
axis, their resulting dynamical evolution during super- 
nova would be approximately in the adiabatic-runaway 
transition region ('I' ~ 0.1 — 1). The result is that their 
pre-Supernova eccentricities (which were probably nearly 
zero due to tidal circularisation) could have been excited 
by a few hundredths to a few tenths, but not by enough 
to have suffered ejection. Although such values fit the ob- 
servations, the system is significantly complicated by the 



mutual interactions amongst all three planets, including 
a resonance locking. Instead, the obs erved pulsar plan- 
ets may be second-generation planets l|Peretg|20ld ). i.e., 
captured (or even formed) after the supernova occurred. 



4.6 Stellar Properties 

Other questions to consider focus on the star itself. How 
does non-constant multi-phase mass loss affect the re- 
sults here? Nonisotropic and/or asymmetric mass loss 
may have a drastic infiu ence on the res ulting cometary 
("Parriott & Alcock"l998') and planetary l|Namounill2005l : 
Namouni & Zhou 2006) evolution. In this case, the sys- 
tem no longer conserves angular momentum, and new 
equations of motion must be derived. How do short 
bursts, periodic or not, of ejected mass accompanying 
pulsating stars affect the planetary orbit? In this case, 
planetary evolution may even undergo several transitions 
between the adiabatic and runaway regimes. The expan- 
sion and/or contraction of the stellar envelope and the 
resulting tidal effects on surviving planets could also play 
an important role in some cases. Tides will compete with 
planetary ejection and possibly eccentricity excitation. 
Further, planets could be expanding their semimajor axes 
- and their Hill Spheres - as they are experiencing tidal 
effects and competing with the expanding stellar enve- 
lope. Some exoplanets will likely be evaporated while 
others will travel through the stellar envelope, accreting 
mass and being subject to a possible non-isotropic mass 
distribution of the stellar envelope. 



5 CONCLUSION 

The variable-mass two-body problem allows for the bod- 
ies to become unbound or highly eccentric. The impli- 
cations of this physical principle affect all dying stellar 
systems which contain any orbiting material. Many Oort 
clouds and wide-orbit planets will have their orbits dis- 
rupted. The extent of the disruption depends crucially 
on their initial semimajor axes, eccentricities, and true 
anomalies, and the subtleties of stellar evolution. Stars 
with progenitor masses of 4Mq — 8Mq will readily eject 
objects that are beyond a few hundred AU distant, and 
excite the eccentricities of the remaining bound mate- 
rial at that distance. Supernovae which produce neutron 
stars eject nearly but not all orbiting material. Con- 
versely, other exotic systems, such as those with black 
holes, could have easily retained planets during their for- 
mation. Stellar mass loss might be the dominant source of 
the free-floating planet population, and orbital properties 
of currently observed disrupted planets in aged systems 
may be tracers of the evolution of their parent stars. 
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